%Minimizing the Variance of Returns
%This m file calculates the weights for the simulated portfolios
% The results of these simulations are reported in Tables 8-11

%The data should be in the following format:
%columns 1-2: blank
%columns 3-8: global fund returns
%columns 9-21: specialized fund returns
%column 22: benchmark index return
%columns 23-36: blank
%column 37: simid ("s" in the code)
%column 38: t (time period)

data1 = textread('data.out', '', 'delimiter', ',','emptyvalue', 1000);

sim=data1(:,37);
[r,c]=size(sim);
dataF = [];

% Loop over the number of simulated portfolios (vary with the data
% frequency and type of simulation):
for s=1:165

%Defining the initial period for the simulations
%Daily data: 50 business days
%Weekly data: 10 weeks 
%Monthly data: 5 months
    period = 50;
    t0=period;                  
    
    t = find(sim == s);
    b = min(t);
    e = max(t);
    
    g= data1(b:e,3:8);
    g2=[];
    i2=1;
    for i2=1:6
        gmc = max(g(:,7-i2));
        if gmc == 1000
        g2= data1(b:e,3:8-i2);
        end
    end    
   
    test=[];
    test2=isequal(test,g2);
    if test2 == 1,
        g2=g;
    end

   glb = mean(g2,2);
   g = g2;
   [rg,cg]=size(g);
   clear g2 i2 gmc 
    
   specialized = data1(b:e,9:21);
    reg2=[];
    i3=1;
    for i3=1:13
        rmc = max(specialized(:,14-i3));
        if rmc == 1000
        reg2= data1(b:e,9:21-i3);
        end
    end    
    
    test=[];
    test2=isequal(test,reg2);
    if test2 == 1,
        reg2=specialized;
    end

    specialized = reg2;
    clear reg2 i3 rmc 

    ret = [g specialized];
    [r,c]=size(ret);

X=[];
x=[];
new=[];
xz=[];
A=[];
bb=[];
Aeq=[];
beq=[];
lb=[];
ub=[];
sigma=[];
retx=[];

n = c; 
nt = n+1; 
T = r; 
gf=cg; 

ew=ones(n,1)*(1/gf);
for h=gf+1:n
    ew(h)=0;
end

xz(1:n,1)=0;
xz(1,1)=0;
X=zeros(1,nt);

i=1;
for t=t0:T
    retx=ret(1:t,1:n);

sigma = cov(retx);    
M=mean(retx);         
A=-M;                 
Aeq=ones(n,1)';       
beq=1;                
f=zeros(n,1);         
ub=ones(n,1);         
lb=zeros(n,1);      
bb=-M*ew;  

[x, favl, exitflag, output, lambda] = quadprog(sigma, f, A, bb, Aeq, beq, lb, ub, xz);
new = [x' i];
X = [X; new];
i=i+1

end

xnew = zeros(period-2,c+1);
xnew= [xnew; X];
xnew2=xnew;
xnew2(:,c+1)=[];
ret2=ret;
ret2(1,:)=[];
last=ret2(1,:);
last(1,:)=0;
ret2=[ret2; last];
retnew=ret2*xnew2';
retnew2=diag(retnew,0);
retnew2=[0; retnew2];
retnew2(r+1,:)=[];
rr=38-(23+c)
z=ones(e-b+1,rr)*1000;

dataF = [dataF; data1(b:e,37:38) data1(b:e,3:8) data1(b:e,9:21) xnew2 z  data1(b:e,1:2) retnew2];

end

dlmwrite('data_VAR.out',dataF, ',');

